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Section 1.0 


SUMMARY 


Investigations into transparent conductive polymers were begun. 
Polypyrrole was electrochemically deposited, but the film char- 
acteristics were poor. A proprietary polymer material supplied 
by Polaroid was evaluated and showed promise as a readily process 
able material. 
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Section 2.0 


INTRODUCTION 


This program will develop the analytical methodology for advanced 
encapsulation designs. From these methods design sensitivities 
will be established for the development of photovoltaic module 
criteria and the definition of needed research tasks. 

The program consists of four phases. In Phase I analytical models 
were developed to perform optical, thermal, electrical and struc- 
tural analyses on candidate encapsulation systems. From these 
analyses several candidate encapsulation systems were selected for 
qualification testing during Phase II. Additionally, during Phase 
II, test specimens of various types will be constructed and tested 
to determine the validity of the analysis methodology developed in 
Phase I . 

During Phase III the following items will be covered: 

1. Correction of identified deficiencies and/or discrepancies 
between analytical models developed during Phase I and rele- 
vant test data obtained during Phase II of the above contract. 

2. Improvement and extension of prediction capability of 
present analytical models. 

3. Generation of encapsulation engineering generalities, princi- 
ples/and, design aids for photovoltaic module design. 
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From these items the sensitivity of module performance to various 
material properties will be determined. This study will enable 
the intelligent direction of research into assessment of module 
life potential by analyzing those materials and their properties 
which through aging would most influence module performance. 

In Phase IV a final optimum design based on knowledge gained in 
Phases I, II and III will be developed and delivered to JPL. 
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Section 3.0 


TECHNICAL DISCUSSION 


3 . 1 Transparent Conductive Polymers 

The equipment for the electrochemical deposition of polymeric 
films has been received and set up. Preliminary experimentation 
was begun using pyrrole in an acetonitrile solution with 0.1 M 
electrolyte. The electrolyte initially used was tetraethylam- 
monium tetraf luroborate . Films were deposited in ITO coated 
glass slides which had been masked with tape to leave a 1.0 x 
1.0 deposition area. 

These intitial experiments yielded poor quality films. Deposi- 
tion of the films began at edges and corners, but did not cover 
the entire area. 

SEM micrographs were taken of polypyrrole films grown on indium 
tin oxide coated glass slides. Figures 1, 2, and 3 show typical 
films. The films start depositing at several nuclei and then 
slowly grow outward. The nuclei grow slightly in area and then 
become thicker, rather than growing in area. Figure 1 shows 
two areas (finger-like areas left of center) which are typical. 
Higher magnifications show a highly textured surface. Chang- 
ing the deposition rate did not improve the character of the 
films. To provide a better nucleation site on the substrate sur 
face, glass slides were prepared with a thin evaporated coating 
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Figure 2. POLYPYRROLE ON ITO 350X 
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Figure 3. POLYPYRROLE ON ITO 3500X 
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of palladium. Figures 4, 5, and 6 show a typical film prepared 
on this substrate. Nucleation sites were much closer but the film 
still had a rough texture. Microscopically the sample showed 
great nonuniformity and none of the samples were suitable for 
characterization . 

Preliminary experiments with a conductive polymer supplied by 
Polaroid were encouraging. The proprietary material was spun 
onto several cells and showed good antireflection properties. 
Electrical properties were not measured. An agreement between 
Spectrolab and Polaroid will be necessary to continue work on 
this material. This is the first conducting polymer which is 
available to be used directly by a spin coating technique. 

3 . 2 Photovoltaic Specialists Conference Participation 

Two papers were accepted for presentation at the 16th IEEE Photo- 
voltaic Specialists Conference in Orlando, Florida. 

3. 3 Wrap-Up of Phase III 

Work has been completed on the additional electrical and struc- 
tural tasks. Two reports are included as appendices to this 
document. This work completes the experimental portion of Phase 
III of this Contract. A final report including Phases I-IV 
will be issued in the fall. 
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Figure 4. POLYPYRROLE ON PALLADIUM 35X 


Figure 5. POLYPYRROLE ON PALLADIUM 350X 
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Figure 6. 


POLYPYRROLE ON PALLADIUM 3 50 OX 



Section 4.0 


CONCLUSIONS AND RECOMMENDATIONS 


There are no conclusions and recommendations for this period. 


Section 5.0 
PLANNED ACTIVITIES 


During the next period Spectrolab will continue investigations 
into the electrochemical deposition of transparent conductive 
materials. The preparation of the final report on Phases I-IV 
will begin. 
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SUMMARY 


A method was developed for calculating the magnitude and location of the maximum 
electric field for a family of solar-cell-like shapes. Simple formulas for use by photovoltaic 
module designers were developed. They provide quantitative information on the effects of 
the cell shape, cell thickness, and pottant thickness on the electrical stress intensification at 
the cell edge. It is recommended that a mathematical investigation be performed to verify 
what the design guidelines for advanced thin-film modules should be. 

A method for calculating the lines of force for three-dimensional electric fields was 
developed and applied to a geometry of interest to the photovoltaic program. 
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1.0 INTRODUCTION 


1.1 BACKGROUND 

This report documents the latest in a series of investigations of the electrical isolation 
design of a photovoltaic module. The objective, as shown in Figure 1-1, is to provide 
guidance on how to design the module to prevent breakdown. The onset of breakdown can 
be expected at the location of the maximum electric field, which occurs at the edge or corner 
of a high-voltage cell. Accordingly, the effects of the cell edge/corner radius, the pottant 
thickness, and the cell thickness on the maximum electric field were investigated. 



Figure 1-1. Electrical analysis objective. 

In the first two phases of the investigation, one-dimensional analyses were performed 
(Ref. 1), and the breakdown voltages of typical module configurations were measured (Ref. 
2). The thrust of the electrical isolation analysis performed in the third phase of the 
encapsulation contract was to develop a method of evaluating the multidimensional effects of 
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the module conductor geometry on the electric field, that is, the electrical stress intensifica- 
tion at edges and corners of the solar cells (Ref. 3). The multidimensional electric field was 
calculated with the finite-element-based NASTRAN Thermal Analyzer, using the analogy 
between thermal and electrostatic fields. To evaluate the maximum electric field accurately, 
the NASTRAN output had to be post-processed; a polynomial was fit to the predicted values 
of potential as a function of distance from the cell surface by the method of least squares, and 
the surface field was calculated from the best-fit polynomial. The accuracy of this method was 
verified by comparison with an exact solution for a geometry similar to that of typical solar 
cells. Two sample geometries were analyzed— a square test coupon used in electrical isolation 
tests performed earlier in this contract and a family of disc-shaped solar cells shown in Figure 
1-2. In the figure, R = cell edge radius of curvature, t c = cell thickness, and t p = pottant 



GROUND PLANE 


Figure 1-2. Disc-shaped solar cell. 


thickness. These sample analyses demonstrated that this finite-element method is a useful 
design tool for evaluating candidate module encapsulation designs. They also showed geo- 
metric limits for which the model incurs numerical difficulties: 1) cells having very sharp 
edges (see Figure l-3a), 2) cells much thinner than the dielectric pottant layer (Figure l-3b), 
and 3) cells much thicker than the dielectric pottant layer (Figure I-3c). 

1.2 PRESENT INVESTIGATION 

The present investigation consisted of the following tasks: 

1. The previous work was extended to the limits for which the aforementioned multidi- 
mensional model encounters numerical difficulties. This analysis is described in 
Section 2. 

2. A method was developed for calculating lines of force of three-dimensional electric 
fields for geometries relevant to photovoltaic modules. This method is described in 
Section 3. 

Conclusions and recommendations are presented in Sections 4 and 5, respectively. 
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a. Sharp-edged cell: t /R>1. 
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b. Thin cell or thick pottant: t p /t c > 1. 
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c. Blunt-edged cell or thin pottant: t p /R<1. 


Figure 1-3, Geometric limits for which the finite- 
element model incurs numerical difficulties. 
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2.0 GENERALIZATION AND EXTENSION OF ELECTRICAL STRESS 
INTENSIFICATION CALCULATIONS 


2.1 APPROACH 

The generalized, extended calculations of the electrical stress intensification were ob- 
tained by an evolutionary process. The key elements were exact solutions for limiting cases, 
physical reasoning, and observing the behavior of the finite-element results when plotted in 
certain ways. 

As in Ref. 3, we evaluate the results in terms of the dimensionless electrical stress 
intensification factor y, defined by 



where 

E,,,, = field for infinite parallel plates 
= V 0 /t„, 

V„ = voltage difference between solar cell and ground plane, and 
E m = maximum electric field. 

(In Ref. 4, y is denoted by t\ and is called the utilization factor. In Ref. 5, y is denoted by t) and 
is called the efficiency factor.) 

The cell shown in Figure 1-2 represents equally well a disc-shaped cell or the edge of a 
rectangular cell. It was found advantageous herein to consider the geometry to be two- 
dimensional planar. This simplified the finite-element analysis and enabled the method of 
conformal mapping to be used. 

2.2 THICK-CELL LIMIT 

The first step in the process was the discovery of an exact solution for the thick-cell limit, 
that is, the geometry shown in Figure 2-1. This solution was obtained by Dreyfus (Ref. 6) 
using conformal mapping. 
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Figure 2-1. Configuration of thick-cell limit (tp/t c — »-0). 


The edge shape for which Dreyfus obtained the solution is not exactly a quarter circle. It 
is, however, very nearly a quarter circle, and the difference is insignificant for our purpose. 

The solution is given in Equations (24) and (27) of Ref. 6. Equation (27) has a typographical 
error: a parameter “a." is typed as "a". In our notation, the correct equations are 


•(tM 


where S is related to t p /R by 


1 (1 + S 2 ) 

. , „ s 1 

arcsinh S — 

(i + s 2 r . 


y s 

t i 

S-arctan S 



JLJil. 

7 r R 


1 

S — arctan S 


+ 


arcsinh S — 


(1 + S ) 


20/2 


(2-3) 


Because of the complexity of Equations (2-2) and (2-3), it is difficult to determine the 
asymptotic behavior as the edge becomes very blunt (t p /R ->• 0) or very sharp (t p /R -> °°). 
Therefore we evaluated the mathematical limits of the terms in these equations and obtained 
the following asymptotic formulas: 



(2-4) 


(2-5) 


The solution is plotted on log-log paper as y versus 1 + (t p /R) in Figure 2-2; in this form, the 
ordinate varies between 0 and 1. (The solution is plotted incorrectly in Refs. 4 and 5.) 







Figure 2-2. Solution for infinitely thick cell. 


The asymptotic formula (2-4) is accurate to within 2 percent for t,,/R < 1. The asymptotic 
formula (2-5) is accurate to within 2 percent for t,,/R > 120. Therefore these simple asymp- 
totic formulas can be used instead of the complicated Equations (2-2) and (2-3) for large 
ranges of the parameter t,,/R. We obtained an even simpler formula by observing that the 
complicated curve in Figure 2-1 is very nearly a straight line for t,,/R > 1; this can be verified 
by placing a transparent straightedge over the curve. We thus found that the exact solution 
can be approximated by 


( 2 - 6 ) 


which is accurate to within 4 percent for 1 < t,,/R < 120. Thus one can use either the exact 
Equations (2-2) and (2-3), the asymptotic formulas (2-4) or (2-5), or the approximate formula 
(2-6), depending on the purpose and the required accuracy. 

The importance to the present analysis of the Dreyfus solution was as follows. It provided 
a basis for extrapolating and checking the finite-element results. It is often difficult to 
determine whether a numerical solution is accurate. We suspected that the finite-element 
results became inaccurate as the cell became thin or thick or as the edge became sharp, but it 
was difficult to evaluate where or by how much the numerical solution became inaccurate. 
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Thus the knowledge that the finite-element results must approach the Dreyfus solution in the 
limit as t p /t c -4 0 was very valuable. The parameter grouping 1 + (t p /R) for plotting the 
Dreyfus solution, which is given in Ref. 4, and the shape of the curve in Figure 2-2 provided 
clues to what parameter groupings and what shapes of master curves would generalize and 
correlate the finite-element results. Furthermore, the Dreyfus solution provides a mathemat- 
ical proof that y -4- 0 as t p /R -> we suspected this result on physical grounds but had 
encountered numerical difficulties in establishing this result with the finite-element method. 
This result means that y -4 0 and E m -4 00 as R -4 0. That is, for a fixed pottant thickness, the 
electric field approaches infinity as the edge of the high-voltage electrode (solar cell) becomes 
sharper; this result is well known to high-voltage equipment designers. It also means that 
y -4 0 as t p -4 00 (holding t p /t c fixed at a very small value). That is, for a fixed solar cell edge 
radius of curvature, the edge behaves increasingly like a perfectly sharp edge as the cell is 
moved further from the ground plane. In this limit, however, the maximum electric field 
decreases as the pottant thickness increases, because y varies as a fractional power of t p ; from 
Equations (2-1) and (2-5), we have 


E m = 


V 0 . Vo 


tpy 


-4 


— 0.9468 (^V 73 
V R / 


1 -4- 0 as t p /R -4 00 , holding R fixed. (2-7) 

2.3 EFFECTS OF CELL EDGE SHAPE 

The effects of the cell edge shape were deduced by physical reasoning, along with an 
examination of the finite-element results and Dreyfus' limiting solution. We considered the 
family of geometries having a fixed value of t p /R. The value of y is then a function of the pa- 
rameter R/L, which can have values between 0 and 0.5. The limit R/t c -4 0 is the thick -cell 
limit for which the Dreyfus solution, described in Section 2.2, applies. The value 0.5 corre- 
sponds to a cell having a semicircular edge. Intermediate values of R/t c correspond to cells 
having an edge with a flat portion. 

By physical reasoning, we concluded that the value of y must be the highest for the 
Dreyfus limit and must decrease as R/t c increases. This behavior occurs because the sun-side 
edge of the cell (the edge away from the ground plane) approaches the shade-side edge 
(where the maximum electric field occurs) as R/t c increases, causing the lines of force to 
crowd together. In other words, the sun-side edge influences the shade-side edge, thus 
increasing the maximum electric field, as R/t c increases. Therefore a plot of y versus l + (t p /t c ), 
holding t p /R fixed, must have the shape shown in Figure 2-3. The finite-element results were 
consistent with this shape. Furthermore, physical reasoning indicated that the slope of such a 


= 0.9468 


R V R 


(J, 

V R 
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1 + — 
l c 

Figure 2-3. Shapes of curves for fixed t p /R. 


curve must approach zero as the Dreyfus limit is approached, because the effect of the sun- 
side edge is significant only for relatively large values of R/t t .. This behavior also was noted 
in the finite-element results. 

This curve provided a check on the finite-element results. The numerical error inherently 
underestimated the maximum electric field and thus overestimated y. Once we realized what 
the shape of the curve in Figure 2-3 is, we were sure that any value of y predicted by the fin- 
ite-element method was wrong if it was higher than the Dreyfus value for the same value of 
tp/R. We discarded these obviously wrong predicted values from the presented data. 

2.4 MASTER CURVE FOR SEMICIRCULAR-EDGED CELL 

Cuddihy (Ref. 7) noted that the finite-element results for the semicircular-edged cell gave 
very nearly a straight line with a slope of —0.5 when plotted on log-log paper in the form y 
versus 1 + (t p /t c ), as shown in Figure 2-4. Thus the following simple formula was obtained: 

y ss (l + ' /2 , for R/t c = 0.5 and 0.25 < t p /t c < 100. (2-8) 

This formula has no physical basis. However, it fits the finite-element calculations very 
closely, and it is consistent with power-law formulas obtained analytically for other geome- 
tries (see Equation (2-5) of this report and Ref. 7). 

2.5 GENERAL MASTER-CURVE FORMULAS 

By combining the results of Sections 2.3 and 2.4, we deduced that the shapes of the curves 
in a log-log plot of y versus 1 + (t p /t c ) must be as shown in Figure 2-5. For fixed values of 
tp/R, each curve has the shape shown previously in Figure 2-3 and approaches the straight 
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Figure 2-4. Shape of curve for semicircular edge. 
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1 Figure 2-5. Shapes of master curves. 

line with a slope of —0.5 shown in Figure 2-4 as R/t,. ->• 0.5. For fixed finite values of t,/R, 
physical reasoning indicates that y 1 as t,,/t,. -»■ 0, because the cell edge appears to be 
infinitely blunt in this limit (see Figure l-3c). 

The finite-element results were consistent with these deduced shapes. In fact, the numeri- 
cal results indicated that, for each fixed finite value of t c /R, the curve becomes parallel to the 
straight line for the semicircular edge for t (J / t t . > 2. That is, 

y cx |l + -^-j for t,,/t t . > 2 and t c /R fixed. (2-9) 


- 24 - 






Pursuing this further, we cross-plotted y versus 1 + (t c /R) on log-log paper. For fixed values 
of t p /t c , the finite-element results turned out to be nearly (to within 10 percent) straight lines 
with a slope of — 1/4, as shown in Figure 2-6. The parameter 


has a minimum value of 3, corresponding to the semicircular-edged cell (t c /R — 2). The 
parameter 


has a minimum value of 1, as desired for a curve fit. Therefore, we deduced that 


y cx 




for t,,/t c > 2. 


( 2 - 10 ) 


Thus we obtained the simple formula 


for 2 < t,,/t,.<100 and R/t, >: 0.025. 



( 2 - 11 ) 



For tp/t c < 2 and R/t c < 0.5, the values of y predicted by the finite-element method agreed 
closely with those calculated by Dreyfus' solution for the thick-cell limit (Equations (2-2) and 
(2-3)). The Dreyfus solution appears to blend with the straight lines having slopes of —0.5. 
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The resultant master curves are plotted in Figure 2-7. The values predicted by the finite- 
element method also are shown on the graph. As was mentioned above, the values are given 
approximately by the Dreyfus solution (Equations (2-2) and (2-3)) for t p /t c < 2 and the curve- 
fit formula (Equation (2-11)) for t p /t c > .2. For example, for t p /t c = 0.5 and R/t c = 0.025, we 
have 



Figure 2-7. Master curves for electrical stress intensification. 


The Dreyfus value is y = 0.37, and this value also is obtained from Figure 2-7. As another 
example, for t p /t c = 20 and R/t c = 0.1, Equation (2-11) yields 




(1 + 20 )' 


- 1\>/4 


= 0.13, 


as is obtained from Figure 2-7. 
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2.6 LOCATION OF MAXIMUM ELECTRIC FIELD 


In addition to the value of the maximum electric field, the location of the maximum field 
is of interest to the designer. The locations determined by the finite-element method are 
shown in Figure 2-8 in terms of the angle 6 defined in the sketch in the figure. 



Figure 2-8. Location of maximum electric field. 


For the thick-cell limit (t p / t c 0), the Dreyfus solution yields 0 = 0. For each value of R/t c , 

the finite-element values of 0 decrease with decreasing t p /t c . The dashed curves, obtained by 
extrapolating the finite-element values to t p /t c =. 0, approach the Dreyfus limit, as they 
should; this is a check on the accuracy of the finite-element results. 

To evaluate the limit of 6 as the cell becomes very thin, the values for the semicircular- 
edged cell were plotted versus t c /t p in Figure 2-9. The value of 0 appears to approach 90 
degrees as the cell shape approaches an infinitesimally thin plane. This result is consistent 
with the exact solution obtained by Maxwell for this limiting shape (Ref. 4). 
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Figure 2-9. Location of maximum electric field on semicircular-edged cell — thin-cell limit. 


2.7 APPLICATION OF RESULTS TO DESIGN OF PHOTOVOLTAIC MODULES 

The results presented in Sections 2.5 and 2.6 provide quantitative information on the 
effects of the module design parameters on the magnitude and location of the maximum 
electric field. 

* Effect of pottant thickness : For a given cell shape, the formulas in Section 2.5 enable the 
reduction in the maximum electric field obtained by increasing the pottant thickness to 
be evaluated. From Equations (2-1) and (2-11), the maximum electric field varies 
approximately, for large values of t p /t c , as 

/ t \ —.1/2 

E m oc (-JL) • (2-12) 

• Effect of cell shape : For a given pottant thickness, these formulas enable the intensifica- 
tion of the electric field produced by making the cell edge very sharp or the cell very 
thin to be estimated. Figure 2-7 shows that very small values of y occur only with very 
sharp edges (as in Figure l-3a) or very thin cells (Figure l-3b). 
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2.8 AREA FOR FURTHER INVESTIGATION 

Because the formula (2-11) was obtained by curve fitting the finite-element results, its 
accuracy for very thin cells is uncertain. Therefore an investigation of the asymptotic 
behavior of y in the thin-cell limit would be useful for designers of thin-film photovoltaic 
modules. The objective would be to evaluate mathematically whether the asymptotic depen- 
dence actually is the simple 1/2 power dependence in Equation (2-11) or another 
relationship. 

A combination of several mathematical techniques appears capable of achieving this 
objective. These tools include conformal mapping (for example. Ref. 8) and perturbation 
methods (Ref. 9). 
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3.0 CALCULATION OF LINES OF FORCE OF THREE-DIMENSIONAL ELECTRIC FIELDS 


This section documents the application of the spline method in the calculation of lines of 
force on planes of symmetry of three-dimensional electric fields. 


3.1 INTRODUCTION 

We present the result of the calculation of the lines of force on planes of symmetry for the 
geometry used in the JPL electrochemical corrosion studies (Ref. 10). To our best knowledge, 
no existing method is available in NASTRAN or elsewhere to perform this task. Therefore we 
developed a method to calculate lines of force, given the NASTRAN output consisting of the 
values of the potential at discrete points. The algorithm we developed uses the method of 
cubic splines. However, extensions were made to account for non-rectangular domains and 
free spline boundary conditions. Section 3.2 gives a brief description of the basic equations. 
The numerical method is outlined in Section 3.3. Finally, the calculations for the electro- 
chemical corrosion test geometry are discussed in Section 3.4. 

3.2 BASIC EQUATIONS 

Some of the following material is a summary of Ref. 11. Let the three dimensional 
function </>(x, y, z) be the electrostatic potential. At any arbitrary point P(x, y, z), let ds be an 
element of the line of force passing through P. (Boldface symbols denote vector quantities.) 
Then ds must be parallel to the electric field intensity E at P. That is, 

ds X E = 0, (3-1) 


where E has Cartesian components (0 X , <p y , (The subscripts x, y, z denote partial deriva- 
tives. For example, denotes d<£/<9x.) This can also be written as 


dx _ dy _ dz 
4>x <Py <P, 


(3-2) 
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(3-3) 


Equations (3-2) for the unknowns x, y, and z have two independent solutions: 


^,(x, y, z) = c„ 
^ 2 (x, y, z) = c 2 . 


where c, and c 2 are constants. Therefore, at every point P, there exist two stream functions. To 
obtain the relationship between the electric field intensity and the stream functions and \p 2 , 
we consider the following. The surfaces = c, and \f/ 2 = c 2 intersect along the line of force L, 
(Figure 3-1). The electric field intensity E is normal to both and V^ 2 on L: 


E • V^, = 0, 
E • V\p 2 = 0. 


(3-4) 



Figure 3-1 . Lines of force and stream functions. 


That is, the electric field intensity E is normal to the plane formed by the vectors V^, and 
V^. Hence there exists a scalar function q(x, y, z) 4= \f/ t and 4= \p 2 such that 


q(x, y, z) E = X V^ 2 . (3-5) 

For a region having uniform permittivity and no charges, 

V • E = 0. (3-6) 

Taking the divergence of (3-5) and using the vector identity, 

V ■ (Vf, X V^ 2 ) = 0, (3-7) 

we have 

E • Vq = 0. (3-8) 
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This implies that 

q = constant, 
and we can set 

q = 1. (3-9) 

Thus we obtain 

E = V^, X Vfa. 

For applications where solutions contain planes of symmetry, the equations for the lines 
of force greatly simplify. If we designate z = c as a plane of symmetry, d<t>/d z must be zero on 
z = c. Equations (3-2) for the lines of force on the plane of symmetry z = c have solutions 

iA.(x, y; z = c) = \ fax, y), (310) 

fa(x, y; z = c) = constant. 

The intersections of the surfaces fa = C! and fa = c 2 are the lines of force. Hence, a set of lines 
of force lies on the plane of symmetry. On the plane of symmetry, the partial derivatives, 
dfa/dx, dyp.y/dy, and d^/Jdz are zero. Substituting (3-10) into (3-5) on the plane of symmetry 
with q = 1 yields the orthogonality conditions 


*r- p(x ' y) i7' 

■af-p (x -y)ir 


where p(x, y) 
condition: 


(-dfaldz) 


(3-11) 

Equation (3-11) is equivalent to the familiar orthogonality 


d\p d<t> d\pd<j) 

— =V\b-V<t> = 0 . 
dx dx dy dy 


(3-12) 


In general, on the plane of symmetry (z = c), the second partial derivative with respect to 
z is nonzero: 


aV 

dr 


+ o. 


(3-13) 


Therefore the two-dimensional Laplacian is, in general, nonzero on the plane of symmetry 
z = c: 


d 2 \ p d~\p 

dx 2 dy 2 


* 0. 


(3-14) 
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This implies that the multiplicative scalar function p(x, y) is not constant on the plane of 
symmetry. 

For strictly two-dimensional problems, p = 1. To see this, the Cauchy-Riemann equations 


<90 _ dip 
dx dy 
d 0 dtp 


(3-15) 


dy dx 

imply that both the electrostatic potential and the stream function satisfy the two-dimension- 
al Laplace equation: 

aV . aV _ n 

dx 2 + dy ~ 

d 2 <P , d 2 <p n 
dx- dy 2 


3.3 NUMERICAL METHOD FOR LINES OF FORCE ON PLANES OF SYMMETRY 

If 0 is known everywhere on the plane of symmetry, equipotentials can easily be 
computed; and the lines of force can similarly be determined if 0 is known on the plane of 
symmetry. From now on, we shall assume we are treating the problem where the plane of 
symmetry is given by z = c. 

However, we need to determine 0 if 0 is known only at discrete points. This can easily be 
done as follows. Given 0, the method of splines is used to obtain the numerical approxima- 
tions of the partial derivatives. Knowing the partial derivatives of 0, the orthogonality 
condition (3-12) determines the partial derivatives of 0 up to the multiplicative function 
p(x, y; z = c). Finally, the spline equations for the two partial derivatives <90/dx and <90/dy 
provide a system of linear equations for the numerical approximation of 0 and the multiplica- 
tive function p. 

Section 3.3.1 gives a brief description of the method of splines for rectangular domains. 
An extension to incorporate free spline boundary conditions is given in Section 3.3.2, and the 
treatment of non-rectangular domains is discussed in Section 3.3.3. Next, the numerical 
algorithm for the calculation of the lines of force from the electrostatic potential is given in 
Section 3.3.4. Two illustrated examples are shown in Section 3.3.5 to verify the numerical 
method. 


3.3.1 Computation of Partial Derivatives of Function 0(x, y) 

Given the function values <p Vi = 0(x ; , yj), i = 1,2, .. . ,N, ) = 1,2, .. . ,M, on a rec- 
tangular mesh (x;, yj], and partial derivatives on the boundary of the rectangle: dip/dy(x„ y,). 
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d(p/d y(x;, y M ), i = 1,2, . . . ,N, and d0/<3x(x,, y ), <30/<3x(x N , y j), j = 1,2, .... ,M, the method of 
cubic splines can be used to obtain the partial derivatives and <j> at every point inside the 
rectangle (see Ref. 12, Chapter 4). The method is efficient, accurate and stable. 


Extensions to Free Spline Boundary Data 


If the partial derivatives <30/<3x and d<t>/dy are not prescribed on the boundary, additional 
conditions must be added to the spline equations. Here we give two alternatives: 

Lagrange polynomial approximation — the partial derivatives on the boundary 
can be approximated by differentiating the Lagrange interpolation polynomial 
and evaluating at the boundary point. An advantage is that it preserves the 
order of the spline method; a disadvantage is its lack of efficiency. 


Artificial boundary condition at a node adjacent to the boundary — to obtain a 
near tridiagonal matrix structure, the second derivative difference approxima- 
tion is matched at the node adjacent to the boundary: 


d0 , . d(f) 

dx ’ dx 1 0(x,) — 20(x 2 ) + 0(X| ). 

— - j-3 


(3-16) 


(For illustration, the above equation is written for the special case of uniform 
mesh, h; the left hand side is the centered difference approximation to the first 
derivative of <90 /dx at x = x>, while the right hand side is the finite difference 
approximation of the second derivative of 0 at x = x>.) An advantage is it is 
very efficient; a disadvantage of the artificial boundary condition is that it is 
only second order accurate. 

Both alternatives require that the number of nodes, N and M, be greater than 3. 
This is not an altogether unreasonable constraint because many practical prob- 
lems use more than 3 nodes for their numerical approximations. 


3.3.3 Extension to Nonrectangular Domain With Quadrilateral Subelements 

The spline method can be used for arbitrary nonrectangular domains by means of a 
change of variables. Specifically, let the physical plane coordinates (x, y) be functions of the 
transformed plane coordinates (r, s). The standard method of splines with free spline 
boundary data extension (Section 3.3.2) can be used to compute the partial derivatives of x, y 
and 0 with respect to r and s in a rectangular domain in the rs plane. The partial derivatives 
of 0 with respect to x and y can be evaluated using the chain rule: 

d0 dx <30 dy _ <30 

dx dr dy dr dr 
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dip 

dx 

dip 


dip 


dx 

ds 

dy 

ds 

ds 


Solving for dtp/dx and dip/dy: 




/ dip 

d y _ 

dip 

dy \ 

dip 


\ dr 

ds 

ds 

dr ) 

dx 



D 





/ dip 

dx 

dip 

dx \ 

dip 


\ ds 

dr 

dr 

~ds) 

dy 



D 




where 

_ dx dy dx dy 

dr ds ds dr 


(3-17) 


3.3.4 Computation of \p . 

Knowing the partial derivatives of ip at every point inside the nonrectangular domain, this 
section describes the final steps to obtain 
The orthogonality conditions (3-11) are 


dip 

dx 


dip 

P(x - y) ip 


d\p dip 

IT" -p (x ' y) ir' 


where the scalar multiplicative function p(x, y) is to be determined. (For the two-dimensional 
Laplace equation, p(x, y) = 1, and the more efficient method of line integration can be used 
to compute the stream function ip. ) 

The partial derivatives of ip with respect to r and s in the transformed plane are given by 


dr~ 

d\p 

ds 


P(X'Y) ( 
P(X'Y) ( 


dip dx _ dip dy 

dy dr dx dr 

dip dx dip dy 

dy ds dx ds 



(3-18) 


For uniqueness of the solution, we observe that, if \p(x, y) is a function which satisfies the 
orthogonality condition (3-12), then ^ = a -F b ^(x, y) also satisfies the orthogonality 
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condition, where a and b are arbitrary constants and b is nonzero: 

d<t> dp , dp dip d<j> d , • , , dp d 

= (a + b w) + (a ■+■ b \I/) 

dx dx dy dy dx dx V dy dy [ 


= }, / ^0 dp + dp_ dp \ 
\ dx dx dy dy / 


For uniqueness, we impose the following initial conditions: 

'Pi r„ s,) = 0, 

p(r„ s,) = a, 


(3-19) 


where a is a nonzero constant. 

The requirement that the second partial derivatives d‘p/dx 2 and d 2 p/d y 2 of the cubic spline 
polynomial approximation of p be continuous at the interior nodes (x ; , yj), i = 2,3, . . . ,(N — 1), 
j = 2,3, . . . ,(M~ 1), provides 2(N — 1) (M — 1) linear equations. Furthermore, the application 
of the artificial boundary conditions (3-16) and the initial conditions (3-19) provide an 
additional 2(N + M— 1) linear equations. Adding, we have 2MN linear equations for the 2MN 
unknowns: p-,j and p ;j , i = 1,2, . . . ,N, and j.= 1, 2, . . . ,M. 

The system of linear equations can be put into the following form: 


^pix = 


(3-20) 


Here A n and A ( .»are block diagonal matrices with its elements being near-tridiagonal matrices 
of the form 


xxx 

xxx 

xxx 


(3-21) 


xxx 
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A 2I and A 22 are near-block-tridiagonal matrices having the same form as (3-21), with its 
elements being diagonal matrices; and 


Pk — (Pk,1/ Pk,2f • • • /pk,N) T 
<lk = OAk.l/ '/'k , 1 ' ■ • • /’/'k.N) 


(3-22) 


for k = 1,2, . . . ,M. The right-hand side f is a (2NM)-vector given by 

i = 1, 
otherwise. 


1 “ 
I 0 


(3-23) 


where a is a constant. 

The solution algorithm of the system of linear equations is straightforward. The matrix 
is factorized into block LU-form: 


A = 




(3-24) 


where 


Lj, A 2) A,, , 

U 22 = A 22 L 2 i A 12 . 


(3-25) 


Since the submatrix A,, is a block diagonal matrix, its submatrices can each be factorized into 
LU-form by using a very efficient near-tridiagonal matrix algorithm. The submatrix L 21 can 
then be very efficiently solved. The near-block-tridiagonal matrix U 22 is further factorized 
into LU-form using a band matrix solver with partial column pivoting. Finally, the solution 
can be obtained in the usual way: 


I 

L 2 | 

\ 

t 

An 

0 


0 

I 


y = f. 


(3-26) 


A,,\ 

uj x = y - 


3.3.5 Numerical Examples to Verify Method 

In this section, we give two examples. The numerical program developed requires the 
following inputs: 

• on the plane of symmetry, z = c, the physical coordinates (x ; , yj), i = 1,2, . . . ,N, 
j = 1,2, . , ,M(N, M greater than 3); 


the potential <t> (x ; , yj). 
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The output from the program is the numerical approximation of ^ at the input physical 
coordinates. 

In each example, the plane of symmetry is z = 0. On this plane, dtpjd z — 0, but d 2 <j>/dz 2 5* 
0, so that 


d~ <t> d 2 <p 

~dx r + ~ty T 


* 0 . 


Example 1: Circular equipotentials 
</> (x,y,z) = x 2 + y 2 - 2z 2 . 

The equipotentials in the plane of symmetry (z = 0) are circles centered at the origin. The 
lines of force are straight lines passing through the origin; this point is therefore a singularity 
of the stream function. The numerical approximation of ^ and the potential in the region 
1 < x <2 , 2< y <3, are given in Figure 3-2. 



Figure 3-2. Equipotentials and lines of force for circular 
potentials of example 1. 
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Example 2: Potential on a curved domain 

<f> (x,y,z) = exp( — x — 0.5y) cos (\/l.25 z) 

for x and y lying in the region given by 

x(r,s) = r(r + 2)/3 + s(s+l)/2 
y(r / s) = 2 — [r(r+3)/4 -s(s + 2)/3] 

for both r and s between 0 and 1. 

This example demonstrated the validity of the theory and the numerical technique given 
above. In Table 3-1, we have the given <f>, the computed ^ and the computed p. In Table 3-2, 
we have the partial derivatives <J> X and <£ y , computed using the numerical method discussed 
above; also the partial derivatives \j/ t and i p y shown in the same table are computed using the 
same numerical procedure with the input obtained from the output of the numerical algo- 
rithm discussed in Section 3.3.4. The rightmost column in Table 3-2 gives the residual 
+ where the partial derivatives are the values obtained from cubic splines. To 
within the accuracy of our numerical method, the computed ^ is orthogonal to the given 0. 
The exact equipotential and the computed lines of force are shown in Figure 3-3. 



Figure 3-3. Equipotentials and lines of force for example 2. 

S> 
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TABLE 3-1 . THE EXACT 4 ; f AND p(x,y) COMPUTED FROM ALGORITHM 
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0. 1 46008+0 1 
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0. 1 OOOOE+Ol 
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0 12000E+00 

0. 2 1 467E+0 1 

0. 30321E+00 

0. 1 1 2888-01 

-0. 75136E+00 

0. 26667E+00 

0. 1 9867E+0 1 

0 28365E+00 

-0 52079E-02 

0 22379E+00 

0 44000E+00 

0 1 80678+01 

O 26097E+00 

0 46182E-01 

0. 1 1708E+01 

O 640008+00 

0. 1 6067E+0 1 

0 236 1 4E+00 

0 1 5078E+00 

0 14225E+01 

0 86667E+00 

0 1 38678+0 1 

0 210148+00 

0 22881E+00 

0. 49600E+00 

0 1 1 2O0E+G 1 

0. 1 1 467E+01 

0 1839 IE +00 

0 1 9692E+00 

-0 16133E+01 

0 28000E+00 

0 23200E+0 1 

0 23693E+00 

0 38292E-01 

-O 1 4458E+0 1 

0 42667E+00 

0 2 1600E+01 

0. 22165E+00 

-0. 61B00E-02 

-0. 20786E+00 

0 60000E+00 

0 198008+01 

0 20393F ♦ 00 

0 20573E-0 1 

0. 1 1 800E+0 1 

0 80000E+00 

0. 1 7800E+01 

0 18452E +OO 

0 1 1 839E+00 

0. 19230E+01 

0 102678+01 

0 1 56008+01 

0. 16420E + 00 

0. 215118+00 

0 1 1447E+01 

0 1 2800E <-01 

0 1 3200E+0 1 

0 1 4370E+00 

0 21952E+00 

-0 13106E+01 

0 480008+00 

0 25200E+O1 

0 17552E+00 

0 77270E— 01 

-0. 24880E+0 1 

0. 626678+00 

0 23600E+01 

0. 1 6420E+00 

0 49248E-02 

-0 10780E+0 1 

O BOOOOE+OO 

0 2 1800E+01 

0 1 5107E+00 

0, 888 1 0E-03 

0 98297E+00 

0 1 OOOOE+O 1 

0 1 9800E+0 1 

0 13670E+00 

0 82985E-0 1 

0 25472E+0 1 

0 12267E+01 

0 1 7600E+01 

0 12164E+00 

0 19242E+00 

0 2 1 560E+0 1 

0. 1 48OOE+0 1 

O 1 52008+01 

0. 10646E+00 

0. 23261E+00 

-0. 69358E+00 

0 72000E+00 

0. 27467E+01 

0. 1 232SE+00 

0 1 2734E+00 

-0 39624E+01 

0. 866678+00 

0 25867E+01 

0. 1 1 533E+00 

0 31001E-01 

-0 26733E+01 

0. 104008+01 

0. 2406 7E+0 1 

0 10610F+00 

-0 859618-02 

0 24581E+00 

0 1 2400E+01 

0 220678+01 

0 96007E-0 1 

0 47698E-01 

0. 32068E+0 1 

0 14667E+C1 

0. 1 9867E+0 1 

0 85435E-01 

0. 16128E+00 

0 371 50E+01 

0. 17200E+01 

0 17467E+01 

0. 74770E-01 

0 234 19E+00 

0 57047E+00 

0 1 OOOOE+Ol 

0. 300008+01 

0 82085E-0 1 

0 1 8 1 96E+00 

-0. 56632E+01 

0. 1 1467E+01 

0 28400E+01 

0 76791E-01 

0 72874E-0 1 

-0. 55612E+01 

O 132008+01 

0 26600E+0 1 

0 70651E-01 

-0 33360E-02 

-0 17128E+01 

0 1 5200E+0 1 

0 2 4600E+O 1 

0 63928E-01 

0 1 6772E-0 1 

0. 36455E+01 

0 17467E+01 

0 22400E+0 1 

0 5688SE-01 

0 1 23 1 7E+00 

0, 61634E+01 

. 0 . 200008+01 

0. 2 OOOOE+Ol 

0. 49787E--0 1 

0 22279E+00 

0. 31 500E+01 





TABLE 3 - 2 . PARTIAL DERIVATIVES COMPUTED USING CUBIC SPLINES 
AND THE RESIDUALS = 


^■ ncijriicnnjfnnn < r « t « roj ^» nnininrii ^ in < i > o > of\i 

^^^^^,^^ 00000 - 00000 - 00000 - 

I ! I I I I I ( I I I I I I I I I I I I I I I I I 
UuiliJUJlliUJUiUiUJtUUJUJlliUlUJUiliiUJUJlilUiUiUJUJUJ 
- 4 jccDniO ‘^' 0 ' OnninajQ 3 Mno - o -^ rijinO‘no 
n « rNOonNCorucv « , n < r - r ) ONnj - oo*-<roo 
mru.nN<r->o(Dnnnoruo*rua)in — njmch^oNn 
CD ' OO -' ON ’-' O ' Ofjm - aj ' ON - mo'-ooNNMnino 
— ru — — — n^om — ru — — — O'ru*' — d - — noonm — 


mm - o > 0 ' 0 ninin > 0 ' 0 >o 
00000-00000 
i I I I I I I I I I I 
UJUJUJUJUJUJliJUJUJUJUJ 
roo<*CD — na-'Ocvinn 
aj (> N ( hfviNa 3 ino ^ — 
nin* 4 j*oocoa)<ro<rao 
n-ooNn-n-ruch 
in ^ru--O'0-n-- 


oooooooooooooooooooooooooooooooooooo 

i i i i i i i i i i i i i i i tit 


0 0-0 
o o o o 

4- 4- ! 4- 

i Id ui UJ UJ 
f) -C «f N 
N pv 00 'O 
N T [MD 
O N O' 

n oj - n 


o-o 

o o o 

+ l +■ 
LU UJ Ui 

n >o n 

ID rf U*) 

n in 
CcJ >o o 
oj >o n 


i o o o 
i o o o 

+ 4- + 

I UJ UJ UJ 

o cu - 
<r n 
CD *0 >0 

cd n n 
— — 


o o o 
o o o > 
♦ + + 
UJ UJ UJ I 
CD >0 O ' 

n o - i 

CD CD fiJ I 

f ^ ■0 I 

— n cvi i 


o-o 

I o o o 

+ I 4- 
I UJ UJ UJ 

n ru o 
o in n 
CO - QD 

O >0 o 

nr«in 


I o o o o 

I o o o o 

4- 4- 4- 4- 

I UJ UJ UJ UJ 

CD O CD (h 

o CM CD m 

o — N ID 

ra n in 
— cu n — 


0 o o o 

1 i t 


O O O 
I I 


0 O o 

1 I I 


o o o < 

I I 


0 o o 

1 I 4- 
UJ UJ UJ 1 

n >o 
o m n 
*n co <r 
ru o n 

0>N h 


0 0-0 
o o o o 

4- 4- I 4- 

I UJ UJ UJ UJ 

i n *0 o N 

i o n in o 

o o ru - 

o — in m 

ru — co — 


- o o 

0 o o < 

1 4- 4- 
UJ UJ UJ ! 
O O O I 
D^ in 

CD O N I 
n ru n 
ru — — i 


o - - 
o o o 

+ i i 

UJ UJ UJ 

>o in O' 
ru o n 
go ru - 

— CD 

ru CD N 


oo-oo-oo—oo—oo— 

ooooooooooooooo 

T+I++I++I++I++I 

UJUJUJUJUJUJUJUJUJUJUJUJUJUJUJ 

no — ni'00'N«r>oaj-coin-n 
Ihofr - Notciioooaj'ouift 
n - ODin <*-< ro ‘<* 0 ‘ Oin «* tu - 
Nn'O^mninin - ru — a* — n n 
— — oru — — — — rururum — — n 


o o o o 

t > 


0 o o 

1 


0 o o 

1 I 


ooooooooooooooo 

III III 


00000000000—000 — — — — — — — — — — — — — — — —— — — — — 
oooooooooooooooooooooooooooooooooooo 

*f*4*f + + + + + + +t+ + +l I I I 1 I II | | I I III II I II I 

UiUJUJUJUiUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUJUIUJUJUJ 
-^■noruN-cDNNNrU'ONru-NGDinchNnininoo^ru^nooincoin- 
oocuoruru«tmrucDcoruin(Do-n-rjnjNChOO>oinnN<on^(><i^NO' 
u^iHo^^rcDru — - ON^ajCDorun — o*NO«raiNrdCDCDCu — (Dinu-NNto*^ 
®Nin«rw— in^n — o — - — ocuru — NramaDon — NtomruNON^r-N* 

— — — O', cd NCOcoN > o > oinoinin < r < j - n < rr r )nncunj 

ooodoooooooooooooooooooooooddooooooo 

i I i I ! I I I I i I I I i l I I I I I I I I I I I I I I I I I I I I I 


O O O O 

coco 


hJ Cu O' O 

m n o 
n n *n n 
- cc m ru 
CD r<j ru CU 


coo 
O O o 

4- 4- 4- 
UJ UJ UJ 

n n 

N O' O' 
o m (h 
*0 n o 
ru ru ru 


o o o 
o o o 

♦ + 4- 
Ui UJ UJ 
N in — 

m ru n 
< r n 

(D o 


o o o o o 

o o o o o 

4- 4- 4- + 4" 

UJ UJ UJ UJ UJ 

— <r - o 

- o -o >o 

in — > o — 

n 'O m n ru 


o o o 
o o o 

♦ + 4- 
UJ UJ UJ 

o cd n 
<r in 
'Onin 

o cu - 


o — — 
o o o 

+ i i 

UJ UJ UJ 

o o n 
ru N CD 
>o — in 
0*0 10 
— O CD 


0 o o o o o o 

1 I I I I I I 
UJ UJ UJ UJ UJ UJ UJ 
-'OOO'ODD- 
CO O* *0 P* — ^ (h 
cd o — o ^ «r o 

— *0 O O *0 O' 
N CD N N *0 in <T 


0000 

1 I I I 


0 O O 

1 I I 


0 o o 

1 I I 


0 o o o o 

1 I I I I 


0 o o 

1 I I 


0 o o 

1 I I 


0 o o o o o o 

1 I I I I I I 


??????????????????<???.??????????????? 
LLiaJUJuJkJUJUJUJUJUJUJUJUJUjUJUJaJUJuJUJUJUJUJUJUJUJUJiiJuJUJUJUJUJUJUJUJ 
OOOOOONNNNNNOOOOOOOOOOOONNNNNNOOOOOO 
OOOOOO'O'O'O'O'O'OOOOOOOOOOOOOO'O'O'O'OOOOOOOO 
c ^- OH )^ o ^ CDOOCD « rru * ocDCD < irurU ' OaDCD ' Oru « j ' CDOO 0 ^ o ^‘' O * o<ro 
C i JD ' 0 ^ njo - 0 ‘ C 0 ' 0 n - n -> Ninninn - 0 ‘ MOMn ^ ni 0 * NOCD * 0 ^ajo 
ru-----ru — — — — — rurd — — - — rururu — — — curururu - — ncvrucucuru 

ooooodoo6o6oooooo6oo66o6ddo666ood6dd 


o o o oo- 
oo o o o o 

■K4- + + 4- + 

UJ Ui UJ UJ UJ Ui 

O N O O N O 

o>coo>oo 
O D o o *0 o 
o ^ ru r j u- o 
o — n in n — 

o 6 6 6 6 6 


o o — — o o o — — — o o — — — — — — — — — — 
oooooooooooooooooooooo 

+ 4 - 4 - 4 , 4 , 4 - 4 , 4 , 4 - 4 - + + + + 4 - 4 ' 4 , + 4 , + 4 , + 
uJ’jJUJUjUJuJuJUJUJUJUJUJUJLUUJUJUIUJUIUJUJUJ 
OONOONOONOONOONOONOONO 
OO'OOO'COO'OOOnOOO-OOO’OOO'OO 

oorucDO'OOorucDO'O^^'Oruo^runj^o 
ooorucorjooru<*ruoocu<rNO — nmNO 

<0 CD — — *0 CD ^ »■< CD — 94 h h ^ — h — ru 

dodo ooddoddodo o 6666666 


4 ] 



3.4 LINES OF FORCE FOR THE ELECTROCHEMICAL CORROSION TEST GEOMETRY 


For the details of the electrochemical corrosion test see Ref. 10. We present here the 
numerical solution \p on the two planes of symmetry indicated in Figure 3-4. The electrodes 
were cylinders having quarter-circular cross sections. The electrode thickness was small 
compared with the electrode radius, and the spacing between the electrodes also was small 
compared with the electrode radius. The planes of symmetry ABCDEFG and AHIJKB are 
perpendicular to each other. The plane ABCDEFG (the z = 0.023 plane) is perpendicular to 
the axes of the cylindrical electrodes and is the midplane of the electrodes. The plane AHIJKB 
(the y = 0 plane) contains the axes of the cylindrical electrodes. The intersection of these 
planes is the line connecting the centers of the electrodes (this line is parallel to the x-axis). 
Point A is the midpoint of this line (V = 1/2), and B is the intersection of this line with the 
high-voltage electrode (V = 1). Point D is in the x = 0 plane at the corner of the high-voltage 
electrode. Point E is in the x = 0 plane at the edge of the dielectric used in the electrochemi- 
cal corrosion tests. Point C is on the high-voltage electrode. Point F is on a line parallel to the 
x-axis passing through E, and G is on a line parallel to the x-axis passing through C. Point K is 
at the edge of the high-voltage electrode, and H is on the V = 1/2 plane on a line parallel to 
the x-axis passing through K. Point I is on the V = 1/2 plane at the edge of the dielectric, and 
J is at the edge of the dielectric on a line passing through B and K. 

Figures 3-5 through 3-7 give the solution for the region ABCDEFG. Figures 3-8 and 3-9 
give the solution for the region AHIJKB. The number shown at a given point in Figures 3-6, 
3-7 and 3-9 is the value of the stream function. The lines of force are the curves connecting 
the points having equal values of the stream function. Note that the y-scale is different from 
the x-scale in Figures 3-6 and 3-7; because of the different scales, perpendicular lines do not 
appear to be perpendicular; thus the lines of force do not appear to be perpendicular to the 
electrode. 

The lines of force — electrostatic potential solution is said to be singular at a point 
Q(x,y;z = c) on the plane of symmetry if = <t> y = 0. At the point Q, the stream function is 
multi-valued. The numerical algorithm presented in Section 3.3 cannot handle this type of 
singularity. Such singularities occur in our geometry at the points D and K in Figure 3-4. We 
do not include the singularities in our computations. 

Because of the acute scaling of the modelled cell, the NASTRAN model has only five 
points on lines parallel to AB in Figure 3-8. Consequently the NASTRAN numerical approxi- 
mation along AB is coarse. Our calculation of the partial derivatives along AB does not 
possess the symmetric behaviour of the exact solution. To compensate for this deficiency, the 
numerical algorithm in Section 3.3.3 is modified: the symmetry boundary condition 



is applied along AB instead of the artificial boundary condition (3-16). 
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Figure 3-4. Electrochemical corrosion test specimen geometry. 
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Figure 3-8. AHIJKB plane — qualitative. 
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Figure 3-9. AHIJKB plane — quantitative. 
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4.0 CONCLUSIONS 


1. Quantitative information on the effect of photovaltaic module design parameters on 
the magnitude and location of the maximum electric field has been developed. 

2. For a family of solar-cell-like shapes, simple formulas for calculating the electrical 
stress intensification at the cell edge have been developed. 

3. The electrical stress intensification is severe only for cells having very sharp edges or 
for very thin cells. 

4. The finite-element method is accurate for modules not having any of the extreme 
geometries mentioned in the preceding conclusion. 

5. For extreme geometries — 1) cells having sharp edges, 2) thin cells, or 3) blunt-edged 
cells or thin pottants — the finite-element results can be extrapolated by means of the 
simple formulas. 

6. A method for calculating the lines of force for three-dimensional electric fields has 
been developed. 
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5.0 RECOMMENDATION 


A mathematical investigation of the asymptotic behavior of electrical stress intensification 
for thin cells should be performed to verify what should be the design guidelines for 
advanced thin-film photovoltaic modules. 
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SUMMARY 


In the first phase of the encapsulation contract, finite element analyses were per- 
formed to predict solar cell stress due to uniform temperature changes of the flat-plate 
terrestrial photovoltaic modules. The results of the analyses were translated into a 
generalized master curve which provided a desktop capability for predicting cell stress for 
any combination of cell, pottant, and structural panel. 

The finite element analyses assumed that out-of-plane bending of the module occurs. 
The degree of out-of-plane bending is unknown. To evaluate the bending assumption, the 
cell stresses were calculated in the present phase of the contract with the assumption that 
no out-of-plane bending occurs. 

For small pottant thicknesses, the cell stress is 10-25 percent higher without bending 
than with bending; this difference is within the overall accuracy of the analysis. For large 
pottant thicknesses, the cell stress is 25 percent lower without bending than with bending, 
for glass and wood structural panels having thicknesses of 0.125 in. and 0.25 in., respec- 
tively. For a 0.08 in. thick steel panel, cell stress is always higher without bending than 
with bending. 

The master curves derived previously for predicting the cell stress can be used as 
presented originally, with the following exceptions: 

• The master curve for solar cell thermal stress analysis should not be used for a cell 
length greater than 4 inches. 

• For beam-like solar cells, the elastic modulus, E, should be input to the master curves 
to compute solar cell stresses. For round or plate-like solar cells, the quantity E / 
(1— v) should be used, where v = Poisson’s ratio. 
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INTRODUCTION 


In the first phase of the encapsulation contract, finite element analyses were per- 
formed to predict solar cell stress due to uniform termperature changes of the flat-plate, 
terrestrial, photovoltaic modules (Reference 1) . The boundary conditions imposed on the 
finite element model permitted out-of-plane bending as well as change in length of the cell, 
pottant, and structural panel laminate. The results of the analyses were translated into a 
generalized master curve which provided a desktop capability for predicting cell stress for 
any combination of cell, pottant, and structural panel (Reference 2) . 

In the current phase of the contract, closed-form analyses were performed to predict 
temperature-induced solar cell stress when out-of-plane bending is not permitted. This 
report presents a comparison of the results of the no-bending closed-form analyses to the 
results of the with-bending finite element analyses in order to evaluate the sensitivity of 
solar cell stress to the out-of-plane bending boundary condition. 

In addition, the accuracy of the solar cell stress master curves for temperature and 
pressure loading of the module are discussed. 

DISCUSSION 

The MSC/NASTRAN finite element model used in the first phase of the encapsula- 
tion program to predict solar cell temperature stress is shown in Figure 1. 

The model consisted of rectangular plate elements which simulated the structural 
panel, the pottant, the cell, and other layers such as front or back covers. As shown in Fig- 
ure 1, only half of a cell was modelled. Symmetric boundary conditions were imposed along 
the plane through the cut edge of the cell, and free-edge conditions were imposed along the 
plane between adjacent cells. In other words, the model behaves as a cantilever beam with 
the left-hand edge considered fixed and non-rotating. 

The boundary conditions imposed on the model permitted the cell and the structural 
panel to interact as a function of their respective thermal stiffnesses , Ea, where E is the 
modulus of elasticity and a is the coefficient of thermal expansion. The interaction 
resulted in bending and stretching of the cross-section, as shown in Figure 2. In the figure, 
shear deformation of the elastomeric pottant is clearly shown. 
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Figure 1. Finite-element structural model for determination of stresses in module 
construction elements in vicinity of a centrally located cell. 



Figure 2. Deflection of structural panel, pottant and cell for 100°C temperature excursion. 



The out-of-plane bending capability was provided to enable deformation consistent 
with that of a classical bimetallic plate when subjected to a uniform temperature change. 
A bimetallic plate consists of two layers of materials which are bonded together and which 
have different coefficients of thermal expansion and elastic moduli. Such a plate bends 
out-of-plane when it undergoes a uniform temperature change (Reference 3) . In the case of 
a photovoltaic module, the two layers are the solar cell and the structural panel, bonded 
together by the pottant. 

The photovoltaic module differs from a classical bimetallic plate since the solar cell 
layer is not continuous, but instead consists of relatively small tiles. As a result, the 
temperature-induced loads in the cell are continuous only within the region of each cell. 
The temperature-induced loads in the structural panel, however, are continuous over the 
entire module surface. It is expected, therefore, that the out-of-plane bending of a photo- 
voltaic module will be less than that of a classical bimetallic plate. Since out-of-plane 
bending was expected to reduce solar cell stress, analysis was undertaken to determine 
solar cell stress when out-of-plane bending is not permitted. 

SOLAR CELL STRESS WITHOUT MODULE BENDING 

A solution for temperature stresses in isotropic laminated plates bonded together by 
flexible adhesives was derived by L.W. Butterworth and R.K. Yasui of JPL (Reference 4) , 
and is presented in Reference 5. The analysis model is shown in Figure 3. 


-I 



Figure 3. Model of three 
laminated layers. 
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It was assumed that the adhesive (pottant) layer had shear stiffness but no axial 
stiffness. The basic equilibrium equations equate the change in axial stress, <r i( in each of 
the outer layers to the shear stress, r, in the bond layer. In the figure x is the distance 
from the center of the joint, and the subscripts 1,2, and 3 refer to the respective layers (in 
our case, the solar cell, structural panel and pottant) . The solutions for axial stress in the 
outer layers, and shear stress in the adhesive, are given by the following equations: 


Ei A a AT 
1 + m n 


cosh 2/Sx ^ 
cosh /3L / 


<r 2 = — m^! 


where 


2Ejt! A a AT , /SL sinh 2/8x 
1 + m n \ cosh /SL 

o' ! = solar cell axial stress 
<r 2 *= structural panel axial stress 
r = pottant shear stress 
ti = solar cell thickness 
t 2 = structural panel thickness 
Ej = solar cell elastic modulus 
E 2 = structural panel elastic modulus 
m = t t /t 2 
n = Ei /E 2 


L = solar cell length 

A a = difference between the coefficients of thermal expansion of 
the solar cell and the structural panel 

AT — change in temperature 

= cell temperature — ambient temperature 

cosh = hyperbolic cosine function 

sinh = hyperbolic sine function 


( 1 ) 

( 2 ) 

(3) 


The dimensionless term /SL is a measure of the relative stiffness of the cell, the 
structural panel, and the pottant, where 
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1 G, 


4 t 3 VEjtj 
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^ 2 t 2 ) 


G 3 = pottant shear stiffness 
t 3 = pottant thickness 


For flexible pottants typical of photovoltaic module construction, 0L < 4. If the pottant is 
very stiff (i.e., high modulus and/or small thickness) , /3L is large. 

Solar cell and pottant stresses for a typical module construction were computed as a 
function of pottant thickness and are compared to the results of the previous with-bending 
analyses in the next section. 


RESULTS 

Solar cell stress was computed for EVA pottant (E = 1000 psi) , a 4 in. X 4 in. X 0.01 
in. silicon solar cell, and glass, wood, and steel structural panels. The results are plotted in 
Figures 4, 5, and 6. 



Figure 4. Cell stress for a glass structural panel, EVA pottant, AT = 100°C, 
4 in. x 4 in. x 0.01 in. silicon cell. 
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Figure 4 is a plot of solar cell stress for a glass structural panel. For pottant 
thicknesses less than 1 mil, the solar cell stress is approximately 25 percent higher when 
bending is not permitted, which was expected. For pottant thicknesses greater than about 
6 mils, the cell stress is higher when bending is permitted. This is probably due to the ef- 
fect of the offset of the solar cell from the neutral axis of the laminate. As the pottant 
thickness increases, the offset bending influence becomes more significant. For a pottant 
thickness of 15 mils, the cell stress is approximately 25 percent higher when bending is 
permitted. 

Similar trends are seen in Figure 5 which is a plot of solar cell stress for a wood 
structural panel. For pottant thicknesses less than 1 mil, the cell stress is about 25 percent 
higher when bending is not permitted. For pottant thicknesses greater than about 3 mils, 
the cell stress is higher when bending is permitted. For a pottant thickness of 15 mils, the 
cell stress is approximately 50 percent higher when bending is permitted. 

Figure 6 is a plot of cell stress for a steel structural panel. In this case, the cell stress is 
always higher when bending is not permitted. For small values of pottant thickness, the 
cell stress is about 10 percent higher when bending is not permitted. The difference 
increases as the pottant thickness increases. For a thickness of 15 mils, the cell stress is 
about 70 percent higher when bending is not permitted. 

A possible cause for the difference between the results for the steel panel and those for 
the glass and wood panels is the thickness of the panels. The steel panel (t = 0.08 in) is 
much thinner than either the glass panel (t — 0.125 in) or the wood panel (t = 0.25 in) . 
The smaller thickness of the steel panel probably reduces the significance of the bending 
behavior. 

Figure 7 is a plot of the pottant shear stress for glass, wood, and steel structural panels 
when bending is not permitted. These stresses were computed using Equation 3. Finite 
element predictions of pottant shear stress are not available. The allowable pottant shear 
stress is 100 psi, shown by the dotted line. It is seen that pottant shear stress is highest 
for steel and lowest for wood. The allowable shear stress of 100 psi is exceeded for pottant 
thicknesses less than 3 mils for a steel panel, and pottant thicknesses less than 2 mils for a 
glass panel. Pottant shear stress is always less then the allowable for a wood panel. It 
should be noted that the pottant shear stress is zero at the center of the cell and highest at 
the edge of the cell. 
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200 


EVA POTTANT AT = 100°C 



ADDITIONAL CONSIDERATIONS 
1. Solar Cell Length 

In equation 1, it can be seen that as /3L is increased, the cell stress approaches an up- 
per bound given by 

EjAaAT 

(T i = 

1 + mn 

On the other hand, in the previously derived master curve for cell temperature stress 
(shown in Figure 8) it is seen that the cell stress increases as Lc 5 , where Lc is the cell 
length. The finite element analyses from which the master curve was generated did not 
consider cells greater than 4 inches in length. Therefore, the master curve should not be 
used for cell lengths greater than 4 inches until additional analyses are performed to verify 
or modify the master curve cell length parameter. It should be noted that both the no- 
bending solution and the with-bending solutions predict that the cell stress decreases as 
cell length is decreased. 
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Figure 8. Master curve for thermal stress analysis. 
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2. Biaxial Material Properties 

The master curves for thermal stress analysis (Figure 8) and pressure stress analysis 
(Figure 9) are generalized for the elastic moduli of the solar cell, the pottant, and the 
structural panel. If the cell is long and narrow (i.e., beam-like), the ordinary elastic 
moduli, Ej, of the layers should be input to the master curves to compute cell stress. 
However, if the cell is a rectangular plate, or is circular, the quantity Ej / (l-Vj) should be 
used, where i>j is Poisson’s ratio of layer ‘i.’ This is due to the effective stiffening of a biax- 
ially stressed plate by the Poisson’s ratio effect. 
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